Dielectric breakdown and avalanches at non-equilibrium metal-insulator transitions 
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Motivated by recent experiments on the finite temperature Mott transition in VO2 films, we 
propose a classical coarse-grained dielectric breakdown model where each degree of freedom rep- 
resents a nanograin which transitions from insulator to metal with increasing temperature and 
voltage at random thresholds due to quenched disorder. We describe the properties of the resulting 
non-equilibrium metal-insulator transition and explain the universal characteristics of the resistance 
jump distribution. We predict that by tuning voltage, another critical point is approached, which 
separates a phase of "bolt" -like avalanches from percolation-like ones. 

PACS numbers: 77.55.-g, 71.30.+h, 45.70.Ht, 64.60.Ht 



Vanadium dioxide (VO2), when heated or strained, 
displays an insulator to metal transition with intriguing 
non-equilibrium collective behavior, portrayed in a re- 
markable series of recent experiments • Strong elec- 
tron correlations drive the microscopies of this metal- 
insulator transition, where a delicate interplay among 
structural, electronic and spin degrees of freedom takes 
place @. However, as we argue in this Letter, the uni- 
versal features of the observed resistance jumps can be 
understood via appropriate generalizations of previously 
studied models of dielectric breakdown 0, 01 • By tun- 
ing two natural control parameters, the applied voltage 
V and the contrast h (the ratio of conductances of the 
insulating and metallic domains), we show that the ex- 
isting experiments are in the small h regime, where a 
crossover, in small samples, between a low-V^ percolating 
phase and a high-y "bolt" phase takes place. As h be- 
comes larger, this crossover evolves to a sharp transition 
with novel critical properties. 

The VO2 films studied in Ref. [lj had a thickness of 
90 nm, widths ranging from 2 /im to 15 /Ltm and lengths 
ranging from 200 nm to 4 /im. X-ray diffraction stud- 
ies of films near criticality revealed that stable insulating 
grains have an average linear size of 20 nm 

with 

the sample put under an external voltage V, multiple 
resistance jumps were observed near the bulk transition 
temperature [1J. The statistics of these jumps revealed 
a power law probability distribution P(AR) ~ AR~ a , 
with an exponent a ~ 2.45. The resistance jump distri- 
bution depended strongly on the magnitude of the exter- 
nal voltage, with the largest jump scaling linearly with 
the voltage. Further, in the presence of external volt- 
age, elongated conducting clusters have been observed 
through X-ray diffraction whereas in the absence 
of voltage, percolation-like isotropic clusters have been 
recorded with near- field infrared spectroscopy [H, EH • 

Even though VCVs transition properties are domi- 
nated by electron correlations, we argue that the ob- 
served collective phenomena can be explained in a 



purely classical way, consistent with experimental obser- 
vations The large length scales of the domains 
(~15-20 nm) and the small electron mean-free path near 
the transition (~0.26 nm) suggest that coherence effects 
are unimportant and electron transport is predominantly 
classical (Ohm 's law) [l2[- This high-temperature tran- 
sition (^340 K) cannot be interpreted as a quantum 
phase transition, since the observed ~1% lattice dis- 
tortion suppresses any electronic coherence (l3j . The 
thermal loading must be considered quasi-static because 
the loading rate of the experiments (< 3 K/min, Q) is 
much slower than the intrinsic dynamics of the domains 
(~10~ 3 s, [loj[). Also, some experiments at high voltage 
show a large event that repeats in space || and time 
over repeated cycles of thermal loading, while others [l[ , 
for smaller voltages do not exhibit this repetition. In 
our model, we consider a quasi-static model of classical 
resistors in two dimensions with deterministic dynam- 
ics, and with classical, quenched disorder, hence lead- 
ing to reproducible avalanche sequences. The strongly- 
correlated quantum and statistical physics underlying the 
Mott transition is absorbed into temperature and volt- 
age dependencies of our domain dynamics, which could 
be estimated by using DMFT methods [HEH- 



Motivated by previous successful studies of strongly 



correlated electronic systems at finite temperatures [16- 
[l8l |. we propose an extended dielectric breakdown net- 
work model of coarse-grained regions transforming from 
insulator to conductor with random critical temperature 
thresholds. We study the resistance jump distribution 
and make predictions about the exponent a. In addition, 
we study the probability distribution of avalanche sizes 
P(s) ~ s~ r , where s is the number of resistors trans- 
formed in a single avalanche burst. We explain the ob- 
served qualitative behavior at different voltages, and pre- 
dict the existence of two distinct regimes: a) a percola- 
tion dominated regime [2cl | where scaling appears only in 
resistance jumps and avalanches are isotropic and small, 
and, b) a bolt dominated regime, where avalanches are 
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FIG. 1: Universal scaling and avalanches in the 
high-contrast, percolation-dominated regime, (a) The 

resistance-temperature curve shows a multiple-step structure 
similar to the experimentally observed one (b) The resis- 
tance jump distribution acquires a universal form, for different 
contrast parameters and voltages for L = 128. The exponent 
a = 2.7 agrees qualitatively with the experimentally observed 
exponent 2.45. Additionally, the distributions show finite-size 
scaling, demonstrating the presence of a nearby critical point. 
Insets: The largest resistance jump is observed to scale lin- 
early with 1/L at fixed V = 0.1 (inset in (b), as observed 
experimentally [l[), and linearly with V at fixed L = 128 (in- 
set in (a)) . 

highly anisotropic, almost line-like. Finally, we make a 
number of experimentally verifiable predictions regarding 
the behavior of the system in the different regimes. 

In our model each link i of the network, labelled by 
a variable Si, is thought of as a microscopic "grain" of 
linear size at least of the order of the dephasing length 1$. 
It can be conducting (Si = +1) with conductance ac, or 
insulating (Si = — 1) with conductance crj = h ac- The 
variable < h < 1 is the inverse contrast between con- 
ducting and insulating regions. We enforce bi-periodic 
boundary conditions on a diamond lattice (a square grid 
rotated by 45°) and subject it to an external voltage 
V per link. Experimental observations show that the 
threshold temperature in the insulator to metal transi- 
tion decreases with voltage 0, 0, [lij]. We account for 



Q-10" 3 



a 










■ 




■ 

■ >. 

a 

B 


- 


II II II II 


3750, h = 
225, h = 
20, h = 

1 


255/256 ° 

15/16 

1/2 


2 

y 


) 


10 1 


S 10 2 


10 3 



(a) 



350 








300 
250 


_16 

^12 
_ _i 

(O 




L = 32 — 

L = 64 — / 

L = 128 — - / 


A 200 


- V 4 






05 

V 150 




-20 -15 -10 ,,-5 
(V-V C )L 1/V 




100 








50 


















10 12 14 16 18 20 22 24 26 28 
V 

(b) 

FIG. 2: Novel universality at the percolation-bolt 
transition, (a) The probability distribution of avalanche 
sizes shows universal scaling near h = 1; In this limit, 
AR ~ S. (b) As shown in the inset, the mean avalanche 
size at the critical voltage diverges as L —¥ oo, suggesting a 
continuous phase transition. 

this in the model by transforming the resistor at link i 
from insulator to metal when the following condition is 
satisfied, 

T >T? — bVi (1) 

where T is the temperature of the sample, Vi is the volt- 
age drop across the i th link, and T? is the random zero- 
voltage critical temperature threshold which models the 
disorder (ljj . Equation[T]is a linear approximation to the 
observed voltage dependence of the critical temperature 
threshold 111 ; the exact functional form should be 
irrelevant for the universal behavior. 

In this model there are two cases which have been 
studied previously: V = and the limit h — > 1. At 
V = 0, resistors are not coupled and transform sequen- 
tially one at a time as in percolation. The resistance 
jump distribution for percolation, originally studied in 
Ref. [jo| . displays a multifractal structure with a power 
law tail at large jumps decreasing with an exponent 
a ~ 2.7 (the power law tail is shown in Fig. Hfb)). 
As h — > 1, the model can be studied by an explicit 
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perturbation expansion in powers of (1 — h)/h [2lj . 
The voltage Vi across the i th link satisfies the recursive 
equation V, = V - - h)/h] £\ + l)V v 

where Tjj are the lattice Green functions with dipo- 
lar form at long distances. Their general form 
for a n-dimensional hyper-cubic lattice is Yij = 
J d n fc/(27r) n sin ifci sin ifcj cos(k.r.y sin 2 §fc/), 

where ki, kj are the directional wave vector components, 
and Yij is the vector from the center of link i to center 
of link j. Taking e = (1 - h)/h « 1 we obtain plf . 

V t - V = -(eF/2)S,Ty(^ + 1) + 0{e 2 ). (2) 

Thus, in the singular limit of h — > 1, the model maps 
to a disordered, long-range, frustrated Ising model. This 
mapping is intriguing, because it maps a dielectric break- 
down model with non-additive multi-body interactions, 
to a dipolar Ising model with additive two-body inter- 
actions. The dipolar interaction in this singular limit is 
shared with a model [22j of interface depinning in mag- 
netic hysteresis, where their finger-like structures resem- 
ble our bolts. 

We perform numerical simulations, where a random 
temperature threshold Tj C , drawn from the standard 
Gaussian distribution, is assigned to each link. The 
simulation starts with every resistor in the insulating 
state (Si = — 1 V i). The voltage at individual nodes is 
found by numerically solving the Kirchoff equations 25 1 . 
At each step the resistor for which the condition T = 
Tf — b\Vi\ (Eq. [1]) is satisfied at the lowest possible value 
of T, is transformed into metal, and voltages are recom- 
puted for the entire network. The process is repeated 
until every resistor in the conducting state (Sj = 1 V i). 

In the experiments of Refs. [l], [2, El on VO2, h is 
small (about 10 -3 ) and the voltage appears to be low 
compared to the disorder threshold. In this limit, for 
large resistance jumps, shown in Fig. QJa), the distribu- 
tion has an exponent a ~ 2.7 which is very similar to the 
experimental findings reported in Ref. [l| . The structure 
of the resistance-temperature curve shown in Fig. [ljb) is 
also similar to ones reported experimentally. The size of 
the largest resistance jump scales linearly with the ap- 
plied voltage, as reported in Ref. jlj. This dependence 
on the applied voltage stems from the non-additive multi- 
body interactions of our model, and cannot be achieved 
by previously suggested bond-percolation type models [l[ 
where the size of the largest resistance jump vanishes in 
the large system size limit. A more explicit signature of 
percolation would be the observation of the multifractal 
scaling [2(| expected at low resistance jumps, possibly 
below the experimental resolution. 

When the contrast is smaller (h > 1/2), we find that 
the insulator to metal transition occurs in avalanches, 
with several bonds transforming simultaneously at the 
same temperature. For fixed h (near 1), avalanches 
and resistance jumps are linearly related (AR ~ s(l — 
h)/(2L 2 ) for the diamond lattice) and both show power 



laws and universal scaling (sizes shown in Fig. [2]). As 
the external voltage V is varied, the avalanche size dis- 
tribution evolves from trivial (at V — 0, where resistors 
transform one by one) to a power law at a critical voltage 
V CT (h), to again trivial (one giant avalanche) at V S> V CI . 
This behavior is suggestive of a continuous phase transi- 
tion; we analyze the probability distribution of our sizes 
P(s) with the scaling form P(s) ~ s~ r $(s/ £ CT , £/£), 
where £ ~ \V — V C r\~ v is the correlation length. The n th 
moment of the avalanche size distribution scales as (s n ) ~ 
L T(i+n-r)^^ v _ v cl )L l / v ). These scaling forms fit the 
data with good accuracy as shown in Fig. [2] Figure [2ta) 
shows the universal size distribution and Fig. [2jb) shows 
the distribution of the mean avalanche size, and a fit 
to the predicted scaling form. From these fits we get 
l/v = 0.25 ± 0.24, a = 0.8 ± 0.4, r = a = 1 ± 0.2. We 
have also studied other disorder distributions (e.g. T? 
taken from a uniform or exponential distribution) and 
explored other analytical methods (e.g. changing the crit- 
ical range in the fits and analyzing the size distribution of 
spatially connected pieces of the avalanches) all of which 
confirm the presence of critical fluctuations. 

The phase transition identified above separates a per- 
colative phase from a 'bolt' phase as shown in Fig. [3] We 
estimate the phase boundary by a mean-field theory that 
becomes exact in the limit of h — > 1, V — > 0. In this limit 
the local voltage concentration are unimportant and the 
interactions are additive. The avalanches can be modeled 
as a branching process - a grain (bond) turning metallic 
induces a long-ranged perturbation in the voltage field, 
which can result in a few more grains turning metallic, 
ad infinitum. The voltage change, A.V, due to a single 



V(l-h) 



metallic bond at a distance r, goes as AV(r) oc 
[32I ]. Let A be the average number of grains that turn 
metallic due the perturbation caused by one grain, then 
A oc V log L(l — h) /(l + h). The mean size of the result- 
ing avalanche is given by 1 + A + A 2 + . . .. Thus, setting 
A = 1 yields a phase boundary between a phase with 
small avalanches (percolative phase, A < 1), and a phase 
with large avalanches (bolt phase, A > 1). 

Figure [3] shows the phase boundary V = 7.26(1 + 
h)/(l — h), where the prcfactor 7.26 is obtained by fitting 
the simulation data. It is difficult to notice the logarith- 
mic drift in the phase boundary due to limited simulation 
size, however, the mean-field analysis suggests that the 
phase boundary is at V = in the limit of L — > 00. 
Even though the voltage per bond, V, goes to zero, the 
externally applied voltage diverges as L/logL. This is 
analogous to fracture where the stress at failure goes to 
zero, and yet the net applied force at failure diverges in 
the limit of large length scales 28j . The mean-field the- 
ory yields a avalanche size exponent of r = 3/2, which 
is different from the numerically observed value (Fig. [5] 
a), possibly due to the effect of fluctuations. Finally, we 
have checked that the mean-field theory can also be col- 
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FIG. 3: Fractal-looking clusters and phase diagram. 

The colors in the insets reflect avalanches; the first spanning 
cluster is shown in black. The phase boundary and the error- 
bars are obtained by treating the critical voltage, V cr , as a 
free parameter in data collapses (see Fig. [2} . The percolation 
fixed point at h = 0, V = is likely unstable under coarse- 
graining, and we anticipate that there will be a crossover to 
the critical point (h — > 1) behavior for very large avalanches 
even for high contrast. 



lapsed by using scaling forms consistent with the scaling 
analysis discussed previously (Fig. [2]b). 

Even though we believe that the phase diagram shown 
in Fig. [3] is qualitatively accurate, there are other pos- 
sible scenarios that cannot be entirely ruled out. It is 
possible that V is finite at the transition, as suggested 
by the scaling analysis. It is also possible that this is an 
avoided critical point [2{|, i.e. large avalanches reflect- 
ing a crossover to the critical point at h — > 1, V — > oo. 
However, the avalanche size distribution displays a scal- 
ing collapse (cf. Fig. 2) and a power law in a large 
range. Also, the behavior is fairly independent of h for 
0.5 < h < 1, rendering a crossover unlikely. 

Our minimal model can be verified experimentally 
in the following ways: a) For high voltages bolt-like 
avalanches should appear, leading to bolt-like conduct- 
ing clusters. This property has already been observed 
in Ref. 0], where elongated clusters appear in the pres- 
ence of finite gate voltage, whereas such anisotropy is 
absent when V = [HI], b) At low contrast (h > 1/2), 
mean resistance jumps and sizes (measured, e.g., using 
multiple ESM images) should diverge only at a critical 
voltage, with power law distributions r = a ~ 1. An 



approach to this regime should be easier in hydrostatic 
pressure-controlled systems like organic materials in the 
re-ET family. 

In conclusion, we presented a novel model of avalanches 
for the metal-insulator transition in VO2, bringing to- 
gether recent experimental findings, and also making con- 



crete experimental predictions as the relevant parameters 
are altered. We have identified a novel continuous tran- 
sition controlled by long-range interactions which could 
be observed in particular classes of materials that have 
evidently smaller contrast, like organic materials under 
hydrostatic pressure [2^, [30] or bulk V2O3 [3ljj . Another 
possibility for achieving low contrast is by tuning hy- 
drostatic pressure, approaching the metal-insulator Ising 
critical point [3] ■ 
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